sep_smooths_year2
Load all required libraries.
library(tidyverse)## -- Attaching packages --------------------------------------- tidyverse 1.3.1 --## v ggplot2 3.3.5 v purrr 0.3.4
## v tibble 3.1.3 v dplyr 1.0.7
## v tidyr 1.1.3 v stringr 1.4.0
## v readr 2.0.0 v forcats 0.5.1## -- Conflicts ------------------------------------------ tidyverse_conflicts() --
## x dplyr::filter() masks stats::filter()
## x dplyr::lag() masks stats::lag()library(plotly)##
## Attaching package: 'plotly'## The following object is masked from 'package:ggplot2':
##
## last_plot## The following object is masked from 'package:stats':
##
## filter## The following object is masked from 'package:graphics':
##
## layoutlibrary(broom)Read in raw data from RDS.
raw_data <- readRDS("./year2.RDS")Make a few small modifications to names and data for visualizations.
final_data <- raw_data %>% mutate(log_copy_per_L = log10(mean_copy_num_L)) %>%
rename(Facility = wrf) %>%
mutate(Facility = recode(Facility,
"NO" = "WRF A",
"MI" = "WRF B",
"CC" = "WRF C"))Seperate the data by gene target to ease layering in the final plot
#make three data layers
only_positives <<- subset(final_data, (!is.na(final_data$Facility)))
only_n1 <- subset(only_positives, target == "N1")
only_n2 <- subset(only_positives, target == "N2")
only_background <<-final_data %>%
select(c(date, cases_cum_clarke, new_cases_clarke, X7_day_ave_clarke)) %>%
group_by(date) %>% summarise_if(is.numeric, mean)
#specify fun colors
background_color <- "#7570B3"
seven_day_ave_color <- "#E6AB02"
marker_colors <- c("N1" = '#1B9E77',"N2" ='#D95F02')#remove facilty C for now
#only_n1 <- only_n1[!(only_n1$Facility == "WRF C"),]
#only_n2 <- only_n2[!(only_n2$Facility == "WRF C"),]
only_n1 <- only_n1[!(only_n1$Facility == "WRF A" & only_n1$date == "2020-11-02"), ]
only_n2 <- only_n2[!(only_n2$Facility == "WRF A" & only_n2$date == "2020-11-02"), ]Build the main plot
#first layer is the background epidemic curve
p1 <- only_background %>%
plotly::plot_ly() %>%
plotly::add_trace(x = ~date, y = ~new_cases_clarke,
type = "bar",
hoverinfo = "text",
text = ~paste('</br> Date: ', date,
'</br> Daily Cases: ', new_cases_clarke),
alpha = 0.5,
name = "Daily Reported Cases",
color = background_color,
colors = background_color,
showlegend = FALSE) %>%
layout(yaxis = list(title = "Clarke County Daily Cases", showline=TRUE)) %>%
layout(legend = list(orientation = "h", x = 0.2, y = -0.3))
#renders the main plot layer two as seven day moving average
p1 <- p1 %>% plotly::add_trace(x = ~date, y = ~X7_day_ave_clarke,
type = "scatter",
mode = "lines",
hoverinfo = "text",
text = ~paste('</br> Date: ', date,
'</br> Seven-Day Moving Average: ', X7_day_ave_clarke),
name = "Seven Day Moving Average Athens",
line = list(color = seven_day_ave_color),
showlegend = FALSE)
#renders the main plot layer three as positive target hits
p2 <- plotly::plot_ly() %>%
plotly::add_trace(x = ~date, y = ~mean_copy_num_L,
type = "scatter",
mode = "markers",
hoverinfo = "text",
text = ~paste('</br> Date: ', date,
'</br> Facility: ', Facility,
'</br> Target: ', target,
'</br> Copies/L: ', round(mean_copy_num_L, digits = 2)),
data = only_n1,
symbol = ~Facility,
marker = list(color = '#1B9E77', size = 8, opacity = 0.65),
showlegend = FALSE) %>%
plotly::add_trace(x = ~date, y = ~mean_copy_num_L,
type = "scatter",
mode = "markers",
hoverinfo = "text",
text = ~paste('</br> Date: ', date,
'</br> Facility: ', Facility,
'</br> Target: ', target,
'</br> Copies/L: ', round(mean_copy_num_L, digits = 2)),
data = only_n2,
symbol = ~Facility,
marker = list(color = '#D95F02', size = 8, opacity = 0.65),
showlegend = FALSE) %>%
layout(yaxis = list(title = "SARS CoV-2 Copies/L",
showline = TRUE,
type = "log",
dtick = 1,
automargin = TRUE)) %>%
layout(legend = list(orientation = "h", x = 0.2, y = -0.3))
#adds the limit of detection dashed line
p2 <- p2 %>% plotly::add_segments(x = as.Date("2021-06-30"),
xend = ~max(date + 10),
y = 3571.429, yend = 3571.429,
opacity = 0.35,
line = list(color = "black", dash = "dash")) %>%
layout(annotations = list(x = as.Date("2021-06-30"), y = 3.8, xref = "x", yref = "y",
text = "Limit of Detection", showarrow = FALSE))
p1 p2Combine the two main plot pieces as a subplot
#seperate n1 and n2 frames by site
#n1
wrf_a_only_n1 <- subset(only_n1, Facility == "WRF A")
wrf_b_only_n1 <- subset(only_n1, Facility == "WRF B")
wrf_c_only_n1 <- subset(only_n1, Facility == "WRF C")
#n2
wrf_a_only_n2 <- subset(only_n2, Facility == "WRF A")
wrf_b_only_n2 <- subset(only_n2, Facility == "WRF B")
wrf_c_only_n2 <- subset(only_n2, Facility == "WRF C")
#rejoin the old data frames then seperate in to averages for each plant.
wrfa_both <- full_join(wrf_a_only_n1, wrf_a_only_n2)%>%
select(c(date, mean_total_copies)) %>%
group_by(date) %>%
summarize_if(is.numeric, mean) %>%
ungroup() %>%
mutate(log_total_copies_both = log10(mean_total_copies))## Joining, by = c("date", "new_cases_clarke", "cases_cum_clarke", "X7_day_ave_clarke", "Facility", "collection_num", "target", "mean_copy_num_uL_rxn", "mean_copy_num_L", "sd_L", "mean_total_copies", "sd_total_copies", "log_copy_per_L")wrfb_both <- full_join(wrf_b_only_n1, wrf_b_only_n2)%>%
select(c(date, mean_total_copies)) %>%
group_by(date) %>%
summarize_if(is.numeric, mean) %>%
ungroup() %>%
mutate(log_total_copies_both = log10(mean_total_copies))## Joining, by = c("date", "new_cases_clarke", "cases_cum_clarke", "X7_day_ave_clarke", "Facility", "collection_num", "target", "mean_copy_num_uL_rxn", "mean_copy_num_L", "sd_L", "mean_total_copies", "sd_total_copies", "log_copy_per_L")wrfc_both <- full_join(wrf_c_only_n1, wrf_c_only_n2)%>%
select(c(date, mean_total_copies)) %>%
group_by(date) %>%
summarize_if(is.numeric, mean) %>%
ungroup() %>%
mutate(log_total_copies_both = log10(mean_total_copies))## Joining, by = c("date", "new_cases_clarke", "cases_cum_clarke", "X7_day_ave_clarke", "Facility", "collection_num", "target", "mean_copy_num_uL_rxn", "mean_copy_num_L", "sd_L", "mean_total_copies", "sd_total_copies", "log_copy_per_L")#get max date
maxdate <- max(wrfa_both$date)
mindate <- min(wrfa_both$date)Build loess smoothing figures figures
This makes the individual plots
#**************************************WRF A PLOT**********************************************
#add trendlines
#extract data from geom_smooth
#both extract
# *********************************span 0.6***********************************
#*****************Must always update the n = TOTAL NUMBER OF DAYS*************************
extract_botha <- ggplot(wrfa_both, aes(x = date, y = log_total_copies_both)) +
stat_smooth(aes(outfit=fit_botha<<-..y..), method = "loess", color = '#1B9E77',
span = 0.3, n = 155)## Warning: Ignoring unknown aesthetics: outfit#look at the fits to align dates and total observations
#both
extract_botha## `geom_smooth()` using formula 'y ~ x'
fit_botha## [1] 11.55339 11.58605 11.61891 11.65197 11.68522 11.71869 11.75236 11.78616
## [9] 11.82006 11.85408 11.88826 11.92265 11.95726 11.99318 12.03098 12.06995
## [17] 12.10936 12.14850 12.18665 12.22311 12.25549 12.28660 12.31932 12.35180
## [25] 12.38393 12.41561 12.44672 12.47715 12.50679 12.53606 12.56507 12.59335
## [33] 12.62038 12.64565 12.66719 12.68454 12.69948 12.71382 12.72935 12.74787
## [41] 12.77117 12.79961 12.83141 12.86517 12.89949 12.93298 12.96422 12.99184
## [49] 13.02273 13.05283 13.07972 13.11119 13.14456 13.17718 13.20637 13.22947
## [57] 13.24381 13.24943 13.24919 13.24473 13.23768 13.22966 13.21814 13.20044
## [65] 13.17799 13.15225 13.12465 13.09662 13.06960 13.04504 13.02735 13.00193
## [73] 12.97505 12.94666 12.91672 12.88518 12.85200 12.81714 12.77779 12.73280
## [81] 12.68455 12.63548 12.58799 12.54450 12.50741 12.47351 12.43876 12.40431
## [89] 12.37130 12.34088 12.31053 12.27817 12.24561 12.21468 12.18718 12.16493
## [97] 12.14974 12.14187 12.13740 12.13913 12.15042 12.16758 12.18694 12.20480
## [105] 12.21749 12.22132 12.21895 12.21585 12.21255 12.20957 12.20744 12.20629
## [113] 12.20563 12.20507 12.20423 12.20272 12.20017 12.19617 12.18520 12.16470
## [121] 12.13881 12.11163 12.08727 12.06986 12.06350 12.07273 12.08550 12.09879
## [129] 12.12121 12.14943 12.18010 12.20988 12.23543 12.25342 12.26555 12.27620
## [137] 12.28601 12.29561 12.30563 12.31671 12.32948 12.34350 12.35788 12.37254
## [145] 12.38742 12.40245 12.41764 12.43304 12.44867 12.46451 12.48060 12.49692
## [153] 12.51348 12.53030 12.54738#assign fits to a vector
both_trenda <- fit_botha
#extract y min and max for each
limits_botha <- ggplot_build(extract_botha)$data## `geom_smooth()` using formula 'y ~ x'limits_botha <- as.data.frame(limits_botha)
both_ymina <- limits_botha$ymin
both_ymaxa <- limits_botha$ymax
#reassign dataframes (just to be safe)
work_botha <- wrfa_both
#fill in missing dates to smooth fits
work_botha <- work_botha %>% complete(date = seq(min(date), max(date), by = "1 day"))
date_vec_botha <- work_botha$date
#create a new smooth dataframe to layer
smooth_frame_botha <- data.frame(date_vec_botha, both_trenda, both_ymina, both_ymaxa)#WRF A
#plot smooth frames
p_wrf_a <- plotly::plot_ly() %>%
plotly::add_lines(x = ~date_vec_botha, y = ~both_trenda,
data = smooth_frame_botha,
hoverinfo = "text",
text = ~paste('</br> Date: ', date_vec_botha,
'</br> Median Log Copies: ', round(both_trenda, digits = 2)),
line = list(color = '#1B9E77', size = 8, opacity = 0.65),
showlegend = FALSE) %>%
layout(xaxis = list(range = c(mindate - 7, maxdate + 7))) %>% #buffer here
plotly::add_ribbons(x ~date_vec_botha, ymin = ~both_ymina, ymax = ~both_ymaxa,
showlegend = FALSE,
opacity = 0.25,
hoverinfo = "text",
text = ~paste('</br> Date: ', date_vec_botha, #leaving in case we want to change
'</br> Max Log Copies: ', round(both_ymaxa, digits = 2),
'</br> Min Log Copies: ', round(both_ymina, digits = 2)),
name = "",
fillcolor = '#1B9E77',
line = list(color = '#1B9E77')) %>%
layout(yaxis = list(title = "Total Log10 SARS CoV-2 Copies",
showline = TRUE,
automargin = TRUE)) %>%
layout(xaxis = list(title = "Date")) %>%
layout(title = "WRF A") %>%
plotly::add_markers(x = ~date, y = ~log_total_copies_both,
data = wrfa_both,
hoverinfo = "text",
showlegend = FALSE,
text = ~paste('</br> Date: ', date,
'</br> Actual Log Copies: ', round(log_total_copies_both, digits = 2)),
marker = list(color = '#1B9E77', size = 6, opacity = 0.65))
p_wrf_asave(p_wrf_a, file = "./site_objects/wrf_a_year2.rda")#**************************************WRF B PLOT**********************************************
#add trendlines
#extract data from geom_smooth
#both extract
# *********************************span 0.6***********************************
#*****************Must always update the n = TOTAL NUMBER OF DAYS*************************
extract_bothb <- ggplot(wrfb_both, aes(x = date, y = log_total_copies_both)) +
stat_smooth(aes(outfit=fit_bothb<<-..y..), method = "loess", color = '#D95F02',
span = 0.3, n = 155)## Warning: Ignoring unknown aesthetics: outfit#look at the fits to align dates and total observations
#both
extract_bothb## `geom_smooth()` using formula 'y ~ x'
fit_bothb## [1] 10.72711 10.81434 10.90000 10.98408 11.06654 11.14735 11.22650 11.30400
## [9] 11.37995 11.45440 11.52744 11.59912 11.66951 11.73748 11.80228 11.86459
## [17] 11.92511 11.98452 12.04349 12.10273 12.15791 12.21275 12.27661 12.34781
## [25] 12.42246 12.49671 12.56667 12.62847 12.67824 12.71904 12.75644 12.79038
## [33] 12.82077 12.84754 12.86471 12.86919 12.86519 12.85691 12.84857 12.84438
## [41] 12.84856 12.85841 12.86840 12.87847 12.88857 12.89864 12.90861 12.91842
## [49] 12.93604 12.95366 12.96499 12.97869 12.99333 13.00754 13.01991 13.02905
## [57] 13.03356 13.03459 13.03411 13.03191 13.02780 13.02157 13.00931 12.98900
## [65] 12.96324 12.93467 12.90591 12.87958 12.85830 12.84470 12.83336 12.82958
## [73] 12.83008 12.83345 12.83827 12.84312 12.84659 12.84726 12.85080 12.86151
## [81] 12.87600 12.89086 12.90267 12.90804 12.90355 12.89337 12.88284 12.87032
## [89] 12.85418 12.83279 12.80022 12.75474 12.70103 12.64375 12.58759 12.53720
## [97] 12.49726 12.45669 12.41770 12.38860 12.36136 12.33616 12.31315 12.29248
## [105] 12.27433 12.25884 12.26067 12.28405 12.31353 12.33367 12.32901 12.30237
## [113] 12.26806 12.22793 12.18386 12.13768 12.09126 12.04645 11.99043 11.91491
## [121] 11.82847 11.73972 11.65727 11.58972 11.54568 11.49412 11.44944 11.43272
## [129] 11.41680 11.40608 11.40497 11.41789 11.44923 11.50342 11.58505 11.68905
## [137] 11.80529 11.92365 12.03403 12.12629 12.19031 12.23950 12.29173 12.34364
## [145] 12.39184 12.43295 12.46827 12.50115 12.53129 12.55839 12.58213 12.60222
## [153] 12.61836 12.63024 12.63757#assign fits to a vector
both_trendb <- fit_bothb
#extract y min and max for each
limits_bothb <- ggplot_build(extract_bothb)$data## `geom_smooth()` using formula 'y ~ x'limits_bothb <- as.data.frame(limits_bothb)
both_yminb <- limits_bothb$ymin
both_ymaxb <- limits_bothb$ymax
#reassign dataframes (just to be safe)
work_bothb <- wrfb_both
#fill in missing dates to smooth fits
work_bothb <- work_bothb %>% complete(date = seq(min(date), max(date), by = "1 day"))
date_vec_bothb <- work_bothb$date
#create a new smooth dataframe to layer
smooth_frame_bothb <- data.frame(date_vec_bothb, both_trendb, both_yminb, both_ymaxb)#WRF B
#plot smooth frames
p_wrf_b <- plotly::plot_ly() %>%
plotly::add_lines(x = ~date_vec_bothb, y = ~both_trendb,
data = smooth_frame_bothb,
hoverinfo = "text",
text = ~paste('</br> Date: ', date_vec_bothb,
'</br> Median Log Copies: ', round(both_trendb, digits = 2)),
line = list(color = '#D95F02', size = 8, opacity = 0.65),
showlegend = FALSE) %>%
layout(xaxis = list(range = c(mindate - 7, maxdate + 7))) %>% #buffer here
plotly::add_ribbons(x ~date_vec_bothb, ymin = ~both_yminb, ymax = ~both_ymaxb,
showlegend = FALSE,
opacity = 0.25,
hoverinfo = "text",
text = ~paste('</br> Date: ', date_vec_bothb, #leaving in case we want to change
'</br> Max Log Copies: ', round(both_ymaxb, digits = 2),
'</br> Min Log Copies: ', round(both_yminb, digits = 2)),
name = "",
fillcolor = '#D95F02',
line = list(color = '#D95F02')) %>%
layout(yaxis = list(title = "Total Log10 SARS CoV-2 Copies",
showline = TRUE,
automargin = TRUE)) %>%
layout(xaxis = list(title = "Date")) %>%
layout(title = "WRF B") %>%
plotly::add_markers(x = ~date, y = ~log_total_copies_both,
data = wrfb_both,
hoverinfo = "text",
showlegend = FALSE,
text = ~paste('</br> Date: ', date,
'</br> Actual Log Copies: ', round(log_total_copies_both, digits = 2)),
marker = list(color = '#D95F02', size = 6, opacity = 0.65))
p_wrf_bsave(p_wrf_b, file = "./site_objects/wrf_b_year2.rda")#**************************************WRF C PLOT********************************************** #add trendlines #extract data from geom_smooth # *********************************span 0.6*********************************** #*****************Must always update the n = TOTAL NUMBER OF DAYS*************************
extract_bothc <- ggplot(wrfc_both, aes(x = date, y = log_total_copies_both)) +
stat_smooth(aes(outfit=fit_bothc<<-..y..), method = "loess", color = '#E7298A',
span = 0.3, n = 155)## Warning: Ignoring unknown aesthetics: outfit#look at the fits to align dates and total observations
#both
extract_bothc## `geom_smooth()` using formula 'y ~ x'
fit_bothc## [1] 10.62188 10.73343 10.84033 10.94256 11.04007 11.13285 11.22086 11.30419
## [9] 11.38294 11.45710 11.52667 11.59164 11.65199 11.70400 11.74557 11.77911
## [17] 11.80707 11.83187 11.85593 11.88169 11.90494 11.92736 11.95127 11.97228
## [25] 11.99167 12.01071 12.03070 12.05289 12.07857 12.10575 12.13167 12.15642
## [33] 12.18009 12.20277 12.22899 12.26152 12.29788 12.33561 12.37220 12.40519
## [41] 12.43208 12.45174 12.46581 12.47601 12.48405 12.49165 12.50054 12.51241
## [49] 12.52054 12.52804 12.54171 12.55697 12.57322 12.58990 12.60641 12.62216
## [57] 12.63659 12.65779 12.68909 12.72269 12.75082 12.76568 12.76741 12.76224
## [65] 12.75157 12.73675 12.71917 12.70020 12.68121 12.66360 12.64198 12.62151
## [73] 12.59658 12.56882 12.53980 12.51113 12.48441 12.46124 12.43836 12.41250
## [81] 12.38524 12.35814 12.33280 12.31077 12.29364 12.28856 12.29604 12.30769
## [89] 12.31512 12.30992 12.29259 12.27016 12.24403 12.21559 12.18624 12.15739
## [97] 12.13044 12.08702 12.04409 12.01507 11.97953 11.94114 11.90358 11.87051
## [105] 11.84560 11.83251 11.83254 11.84162 11.85545 11.86970 11.88004 11.89711
## [113] 11.93023 11.97309 12.01941 12.06286 12.09716 12.11600 12.11988 12.11497
## [121] 12.10352 12.08774 12.06988 12.05216 12.03682 12.01304 11.98806 11.96948
## [129] 11.94658 11.92115 11.89500 11.86995 11.84780 11.83036 11.81330 11.79327
## [137] 11.77357 11.75749 11.74832 11.74937 11.76391 11.78756 11.81474 11.84693
## [145] 11.88564 11.93236 11.98652 12.04663 12.11278 12.18508 12.26360 12.34845
## [153] 12.43972 12.53749 12.64187#assign fits to a vector
both_trendc <- fit_bothc
#extract y min and max for each
limits_bothc <- ggplot_build(extract_bothc)$data## `geom_smooth()` using formula 'y ~ x'limits_bothc <- as.data.frame(limits_bothc)
both_yminc <- limits_bothc$ymin
both_ymaxc <- limits_bothc$ymax
#reassign dataframes (just to be safe)
work_bothc <- wrfc_both
#fill in missing dates to smooth fits
work_bothc <- work_bothc %>% complete(date = seq(min(date), max(date), by = "1 day"))
date_vec_bothc <- work_bothc$date
#create a new smooth dataframe to layer
smooth_frame_bothc <- data.frame(date_vec_bothc, both_trendc, both_yminc, both_ymaxc)#WRF C
#plot smooth frames
p_wrf_c <- plotly::plot_ly() %>%
plotly::add_lines(x = ~date_vec_bothc, y = ~both_trendc,
data = smooth_frame_bothc,
hoverinfo = "text",
text = ~paste('</br> Date: ', date_vec_bothc,
'</br> Median Log Copies: ', round(both_trendc, digits = 2)),
line = list(color = '#E7298A', size = 8, opacity = 0.65),
showlegend = FALSE) %>%
layout(xaxis = list(range = c(mindate - 7, maxdate + 7))) %>% #buffer here
plotly::add_ribbons(x ~date_vec_bothc, ymin = ~both_yminc, ymax = ~both_ymaxc,
showlegend = FALSE,
opacity = 0.25,
hoverinfo = "text",
text = ~paste('</br> Date: ', date_vec_bothc, #leaving in case we want to change
'</br> Max Log Copies: ', round(both_ymaxc, digits = 2),
'</br> Min Log Copies: ', round(both_yminc, digits = 2)),
name = "",
fillcolor = '#E7298A',
line = list(color = '#E7298A')) %>%
layout(yaxis = list(title = "Total Log10 SARS CoV-2 Copies",
showline = TRUE,
automargin = TRUE)) %>%
layout(xaxis = list(title = "Date")) %>%
layout(title = "WRF C") %>%
plotly::add_markers(x = ~date, y = ~log_total_copies_both,
data = wrfc_both,
hoverinfo = "text",
showlegend = FALSE,
text = ~paste('</br> Date: ', date,
'</br> Actual Log Copies: ', round(log_total_copies_both, digits = 2)),
marker = list(color = '#E7298A', size = 6, opacity = 0.65))
p_wrf_csave(p_wrf_c, file = "./site_objects/wrf_c_year2.rda")keeping in case
#save(wrfa_both, file = "./plotly_objs/wrfa_both.rda")
#save(wrfb_both, file = "./plotly_objs/wrfb_both.rda")
#save(wrfc_both, file = "./plotly_objs/wrfc_both.rda")#save(date_vec_botha, file = "./plotly_objs/date_vec_botha.rda")
#save(date_vec_bothb, file = "./plotly_objs/date_vec_bothb.rda")
#save(date_vec_bothc, file = "./plotly_objs/date_vec_bothc.rda")#save(both_ymina, file = "./plotly_objs/both_ymina.rda")
#save(both_ymaxa, file = "./plotly_objs/both_ymaxa.rda")
#save(both_yminb, file = "./plotly_objs/both_yminb.rda")
#save(both_ymaxb, file = "./plotly_objs/both_ymaxb.rda")
#save(both_yminc, file = "./plotly_objs/both_yminc.rda")
#save(both_ymaxc, file = "./plotly_objs/both_ymaxc.rda")